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Abstract 

Consider the three-dimensional flow of a viscous Newtonian fluid 
upon an abitrarily curved substrate when the fluid film is thin as oc- 
curs in many draining, coating and biological flows. We derive a model 
of the dynamics of the film, the model being expressed in terms of the 
film thickness and the curvature tensor of the substrate. The model 
accurately includes the effects of the curvature of the substrate, via 
a physical multiple-scale approach, and gravity and inertia, via more 
rigorous centre manifold techniques. Numerical simulations exhibit 
some generic features of the dynamics of such thin fluid films on sub- 
strates with complex curvature. 



Contents 



1 Introduction 



1.1 Equations of motion| 4 



1.2 Boundary conditions] 5 



2 The adopted coordinate system 






2.1 The curvilinear coordinate system 




2.2 Free-surface geometry 





*Dept of Mechanical Engineering, University of Delaware, USA. E-mail: 
royOme . udel . edu 

^Dept of Mathematics & Computing, University of Southern Queensland, Toowoomba, 
Queensland 4352, Australia. E-mail: aroberts@usq.edu.au 

*Dept of Mathematics & Computing, University of Southern Queensland, Toowoomba, 
Queensland 4352, Australia. E-mail: simpsonm@usq.edu.au 



1 



11: Introduction 



3 Lubrication model driven by surface tension 



9 

3.1 Re-scale the problem] 10 



372 Perturbation solution! 12 



3.3 Conservation of mass 13 



4 Centre manifold analysis: Gravity and inertia 




4.1 Basis ol the centre manilold . 




4.2 The general lubrication model 





15 



5 Example film flows 



21 

5.1 Qualitative effect of substrate curvature! 21 

23 



572 Corner flow 



573 flow around a torus! 27 



Kefe 



rences 



28 



1 Introduction 



The importance of thin-film fluid flows in countless industrial and natural 
processes has lead to the development of a variety of mathematical mod- 
els and numerical simulations. This has increased the understanding of the 
various complex processes at play, such as the effects of generalized New- 
tonian and shear-thinning rheology of the liquid p7|, of geometric com- 



plexity of the substrate |pS 
cesses 



of multicomponent mixtures and drying pro- 
of advecting 



of surface contamination |30], |3l| and roughness p2 



and diffusing contaminants, and of moving substrates. See also for example 
II [16, [17|, [735], |73j|. Examples of industrial applications include the 



coating processes of autobodies, beverage containers, sheet goods and films, 
decorative coating, gravure roll coating, etc. Physical applications can be 
found in the biomedical field such as the liquid films covering the cornea of 
the eye or protecting the linings of the lungs ||. 

Herein we consider the slow motion of a thin liquid layer over an arbi- 
trarily curved solid substrate. The fluid is assumed to be incompressible 
and Newtonian, constituted of a single component and uncontaminated by 
surface-active material. The substrate is stationary. The effects of substrate 
curvature on the flow of thin liquid layers driven by surface tension were 
first modeled by Schwartz & Weidner P9|| for two-dimensional geometries. 
Short-wavelength irregularities of the fluid's free surface are quickly leveled 
by surface tension forces, but the long-term evolution of the flow is primarily 
determined by the substrate's curvature. Their calculations confirm qualita- 
tive observations of the thinning of coating layers at outside corners and of 
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the thickening at inside corners. In a subsequent work P0| , they studied the 
joint effect of substrate curvature and the presence of surfactant and showed 
that corner defects may be significantly altered by the presence of Marangoni 
forces. 

The model derived herein exploits the thinness of the fluid layer and leads 
to a "lubrication" model of the dynamics whereby the unknown fluid fields 
are parameterised by only the fluid layer's thickness. Thus, a considerable 
reduction of the dimensionality of the problem is achieved when compared 
to the full Navier-Stokes equations. Such thin-geometry or long-wave models 
have previously been derived in many physical contexts, such as in the me- 
chanics of beams, plates and shells |$4], |22| , coating flows [13], ^7], || |[ |35| , 



shallow- water waves viscous fluid sheets |3"o| , etc. On a flat substrate 
the usual model for the evolution of a viscous film's thickness rj, driven only 
by surface tension, is given by the following non-dimensional equation 



dt 3 



where k ~ V 2 r? is the curvature of the free surface of the film. Based upon 
the conservation of fluid and the Navier-Stokes equations, outlined in §§fTI} 
we derive in §0 the following model for the evolution of a film on a curved 



-|V • [t/ 2 (Vk - ±r/ 4 (/d - K) • Vk| , (1) 



substrate: 

d( 
dt ' 

where C — V ~ \^if + \kik 2 rf is proportional to the amount of fluid locally 
above the substrate; K is the curvature tensor of the substrate; hi, k 2 , and 
k = k\ + k 2 are the principal curvatures and the mean curvature of the sub- 
strate, respectively; and the V-operator is expressed in a coordinate system 
of the substrate, as introduced in §Q. This model systematically accounts for 
the curvature of the substrate and that of the surface of the film. It improves 



the model derived by Schwartz & Weidner |29|] and extends it to flows where 
the substrate curvature has a larger effect on the fluid dynamics. 

Modifications to this model are derived, in when gravity and/or inertia 
are significant influences. The analysis is performed using computer algebra 
and is based on centre manifold theory [^] to assure us that all the dynamical 
effects are incorporated into the model. In centre manifold theory competing 
small effects need not appear at leading order in the analysis; thus we obtain 
the flexibility to adapt the model to a variety of parameter regimes without 
redoing the whole analysis. For further discussion of this and other aspects 
of the application of centre manifold theory to low- dimensional modelling see 
the review in f26|| . Here in particular, the various regimes of gravitational 



forcing around a curved surface are all encompassed within the one model, 
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namely equation (pQ|) . The terms affected by fluid inertia appear very small 
in generic situations, see fl52|), and perhaps may be best used to indicate the 
error in the lubrication models for moderate Reynolds numbers. 

In the last section, §[|, we conclude with some numerical simulations of 
flows on curved surfaces, with attention given to the quantitative differences 
found between our model and that obtained in pi . 



1.1 Equations of motion 

We solve the Navier-Stokes equations for an incompressible fluid of density 
p and viscosity p moving with velocity field u and pressure field p. The 
flow is primarily driven by pressure gradients along the substrate and caused 
by capillary forces characterized by surface tension a and varying due to 
variations of the curvature of the free surface of the fluid. However, there 
may also be a gravitational body force, g, of magnitude g in the direction 
of the unit vector g. Suppose the film has characteristic thickness H. We 
non-dimensionalise the equations by scaling variables with respect to: the 
reference length H; the reference time pH/a; the reference velocity U — a/ fi, 
and the reference pressure a/H. Thus, in this non-dimensionalisation we take 
the view of a microscopic creature of a size comparable to the thickness of the 
fluid; later we require that both the substrate and the free surface curve only 
gently when viewed on this microscale. The non-dimensional fluid equations 
are then 

Vu = 0, (2) 
= -Vp+V 2 u + 6g, (3) 

where 1Z = apH/p 2 is effectively a Reynolds number characterising the im- 
portance of the inertial terms — it may be written as UH/v for the above 
reference velocity — and b = pgH 2 /a is a Bond number characterising the 
importance of the gravitational body force when compared with surface ten- 
sion. 

In §^ we assume that the regime of the flow is characterized by a very 
small value of the Reynolds number so that the inertia term TZDu/ Dt is 
neglected in comparison to the viscous forces in the fluid: 

TZ = UH/v < 1 . 

This is the "creeping flow" assumption of lubrication. Later, in §[|, we rein- 
state inertia into the analysis and determine its leading order effects on the 
dynamics. 



K 



du 

sF + uVu 
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1.2 Boundary conditions 

Preliminary statements of the boundary conditions follows. 

1. The fluid immediately in contact with the substrate does not slip along 
the stationary substrate S, that is 

u = on S . (4) 

2. The kinematic boundary condition at the free surface of the fluid states 
that fluid particles must follow the free surface. 

3. The free surface, denoted by T and assumed free of contamination, 
must have zero-shear (tangential stress), namely 

t • ti = r ■ t 2 = on T , (5) 

where f is the deviatoric stress acting across J 7 , and t a are indepen- 
dent tangent vectors to T . This condition assumes a light and inviscid 
medium (such as air) above the fluid layer — a condition which guar- 
antees the continuity of the tangential stress as one passes across the 
interface T . 

4. The normal surface stress at T must account for the surface tension 
whose effect is to create a discontinuity in the normal stress propor- 
tional to the mean curvature of T: in non-dimensional form 

p = — k + f • n onjF, (6) 

where p is the fluid pressure relative to the assumed zero pressure of 
the medium above, and n is the unit normal to T . 

We do not discuss boundary conditions at the lateral extremes of the sub- 
strate as the substrate is assumed to be so large in extent, when compared to 
the fluid thickness, that the dynamics of the fluid film are largely unaffected 
by the edge boundary conditions. A rational method for deriving boundary 
conditions for nonlinear dynamical models such as (0) is explained in ||21|| . 

2 The adopted coordinate system 

The flow is best analysed and solved in a coordinate system that naturally 
fits the curving substrate. Based upon the principal directions of curvature 
of the substrate, we construct an orthogonal curvilinear coordinate system 
in the neighbourhood of the substrate. This natural coordinate system has 
some remarkable properties which make a systematic analysis tractable. 
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Figure 1: The substrate S is parameterised by variables x\ and x 2 - Together 
with the normal distance y, these form an orthogonal curvilinear coordinate 
system in space with unit vectors ei, e 2 and e 3 . 

2.1 The curvilinear coordinate system 

The curvilinear coordinate system is taken to be an extension into space of a 
natural coordinate system of the substrate. We choose, without loss of gen- 
erality, a curvilinear parameterisation (xi,X2) of S for which the parameter 
curves X\ = Constant and x 2 = Constant generate the lines of curvature of 
S ||33|| , shown schematically in Figure [I| In this case, the curvature tensor 
becomes diagonal everywhere, with diagonal components k a .^ 

Conversely, we may view the surface S as being entirely specified by 
its metric coefficients m a (ii,i 2 ) and its principal curvatures k a (xi,x 2 ) as 
functions of the coordinate variables (xi,x 2 ) assumed to generate the lines 
of curvature as the parameter curves of S. 

We also define the triad of unit orthogonal vectors (e 1; e 2 , e 3 ) at all points 
of S: e a is tangent to curves of constant x a >] and e 3 is normal to S. The 
metric coefficients, m a , allow measurement of lengths on S: the arclength of 

1 Latin indices (such as span the numbers 1, 2, 3 and are attached to spatial quan- 
tities, while Greek indices (such as a or (3) span the numbers 1, 2 and are attached to 
substrate or surface quantities. A dash on Greek indices denotes the complementary 
value, that is, a' = 3 — a. Subscripts following a comma indicate differentiation with re- 
spect to the corresponding coordinate. Unless otherwise specified, we do not use Einstein's 
summation convention for repeated indices. 
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curves on the substrate are found from 

ds 2 = (midxi) 2 + {m 2 dx 2 ) 2 ■ (7) 

Note that the unit vectors vary along S and that their derivatives dei/dx a 
are needed to express the equations of motion and boundary conditions in 
the curvilinear coordinates. 

As shown in Figure [l], the two-dimensional substrate S is the locus of 
the endpoints of the position vector rop = X(xi, x 2 ) , for some domain T> of 
(xi,x 2 ). We prescribe a third coordinate, denoted by either x 3 or y, as the 
distance measured along the normal e 3 from a given spatial point P to the 
surface S. Thus the position vector of points P of the fluid may be written 
as 

r = X(x lt x 2 ) + ye 3 (x 1 ,x 2 ) , (8) 

where the endpoint of vector X belongs to the surface S. A large number of 
important simplifications occur in using this particular orthogonal curvilinear 
coordinate system, (xi,x 2 ,y), that naturally fits the substrate. A definite 
example for the flow on a torus is given in §§ |5.3| . 

We denote by rj(t, x\, x 2 ) the fluid layer's thickness at time t and location 
(xi,x 2 ). The free surface T of the fluid is thus represented by the equation 
y = r](t,xi,x 2 ), and the fluid fills the domain < y < r](t,xi,x 2 ). Note 
that in general, sufficiently far from the substrate for example, the chosen 
coordinate system does not lead to a one-to-one mapping between coordinate 
space and physical space. We require that proper conditions between r], k\ 
and k 2 are satisfied so that intersections of normals of S do not occur within 
the body of the fluid. 

The coordinate unit vectors are independent of y and hence are the same 
unit vectors, (ei, e 2 , e 3 ), as defined on the substrate but now defined in space 
neighbouring S. The corresponding spatial metric coefficients are 

h a = m a (l - k a y) , h 3 = 1 . (9) 

By relating the base unit vectors of the rectangular coordinate system to 



those of the curvilinear coordinate system, see [0, p598] or |19[, the spatial 
derivatives of the curvilinear unit vectors are 

e *,i = ~r h hi - ~r h ^ ■ (10) 

"t k=l h 

However, in the special coordinate system used here, all such expressions 
simplify to the equivalent expressions on the substrate and independent of y: 



G Q ' "T K a Tn a & 3 , G« a' 6 q ' 



e 3 ,a = -k a m a e a , e ii3 = . (11) 
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This choice of curvilinear coordinates is needed only for the derivation of a 
reduced-order approximation of the dynamics of coating flows on S. As seen 
in ([!]) for example, the model will be ultimately expressed in coordinate-free 
terms. 



2.2 Free-surface geometry 

The shape of the free surface is critical in thin film flows: fluid surface cur- 
vature variations create pressure gradients to drive the fluid flow. To denote 
a quantity evaluated on the fluid's free surface we generally use an over-tilde 
(as already done in §§ |1.2|) . 

The position of points P on the fluid's free surface T is given by 

r p = X(xi, x 2 ) = X(xi, x 2 ) + r}(t, xt, x 2 )e 3 (x 1 ,x 2 ) . (12) 

Hence, the surface T is naturally parameterised by (xi, x 2 ), and its tangential 
vectors, ti and t 2 , and unit normal vector, n, are given by 

t Q = -q^- = h a e a + 77 iQ e 3 , (13) 
_ = ^ X ~ 2 oc -h2%iei - ~hir] t2 e 2 + hji^ , (14) 

|tx X i2| 

where h a = m a (l — k a rj) are the metric coefficients at the free surface, and 
where r] i0l = drj/dx a . 

At the free surface, y = rj, the kinematic boundary condition must be 
imposed, namely that fluid particles on the free surface remain on it. This 
leads to 

drj ui dr] u 2 dr] 

m =v -h^-h^_ onT - (15) 

To account for surface tension effects, we need to compute the free surface 
mean curvature k in terms of the substrate principal curvatures and the film 
thickness rj. A tractable route is to recognise that the effect of surface tension 
arises through the energy stored in the free surface. Thus its contribution to 
the dynamical equations, through the curvature k, arises from the variation 
of the surface area with respect to changes in the free-surface shape y = 
r)(t,Xx,x 2 ). The free-surface area is A — j dA = J Adx\dx 2 where 



A = Jhjhl + h\r}\ + h\r] 



2 

2 ! 



is proportional to the free-surface area above a patch m\dx\ x m 2 dx 2 of the 
substrate. The effect of curvature of the free surface, k, is determined from 
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the variation of A with respect to r]: 

d (OA 



hih,2K 



SA 

5rj dx\ \dn\ 



d ( OA 



dx 2 \dr)o 



OA 

dr] 



That is, 



K 



1 

+ A 



An approximation is 



V 2 rj + 



dx\ 
ki 



A 



+ 



2 \ m 2 k 2 



dx 2 V A 
hi + r)\ 



m\k\ 



1 — k\T\ 1 — k2T] 



n 



(16) 



(17) 



where it is sufficiently accurate to use the standard form of the Laplacian in 
the curvilinear coordinates of the substrate, 



V 2 r] 



m\ni2 



d ( irt2 drj\ d ( mi drj 



dxi \ mi dx\ 



+ 



8X2 I m 2 9X2 



The approximation (|T7|) arises directly from the form of the variational ex- 
pression and accounts for changes in the free-surface curvature due to the 
finite depth of the film, and variations in the film thickness. 

Note that throughout this paper, V has two different meanings depending 
upon the context of whether it is applied to three-dimensional spatial fields, 
such as u and p, or to two-dimensional substrate fields such as t] and k. 



3 Lubrication model driven by surface ten- 
sion 

We now derive a model for the flow dynamics when the fluid film and the 
substrate vary on a large scale relative to the thickness of the film: because 
of the two vastly different space scales it may be viewed as a multiple-scale 
analysis. In this case, viscous dissipation acts quickly across the fluid to damp 
out all except the slow dynamics associated with such large scale motion. 
With the assumptions that inertia and body forces can be neglected, by 
considering conservation of mass we derive an equation for the slow evolution 
of the film thickness 77. 
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3.1 Re-scale the problem 

We need to re-scale the non-dimensional governing equations. There are two 
characteristic lengths of the problem: a reference length L measured along 
the substrate and a reference thickness H of the fluid layer covering S. The 
length-scale L is thought of as either the scale of the radius of curvature 
of the substrate, or as the scale on which the film thickness varies. We 
generally expect both to be of the same order of magnitude as they both 
are inversely proportional to substrate gradients, Vn and respectively. 
We denote the ratio H/L by e, and assume e < 1 to be consistent with 
the thin-film/large-substrate assumption. In particular, we consider the flow 
to be non-dimensionally of thickness 1, and so the substrate scale is non- 
dimensionally of large size 1/e. 

With the above in mind, we introduce the scaled curvatures (temporarily 
indicated by a "*" superscript) 

k a = e/c*, and k = ek*, (18) 

to express that the substrate and free-surface curvatures are O (e). Then we 
scale substrate coordinates and metric coefficients according to 

m a = ^-to* , x* a =x a . (19) 

This form is appropriate if the substrate coordinates are naturally non-dimen- 
sional, such as the angular latitude and longitude coordinates on a sphere. 
(If you consider the substrate coordinates as naturally lengths, then the al- 
ternative scaling m a = to* , x* a = ex a is appropriate; in this case you would 
consider slow-space scale.) The re-scaled spatial metric coefficients 

then become 

h* a = m* a (l-€k* a y). (20) 

Seeking a slow flow leads to the following re-scaling of the nondimensional 
variables for pressure, velocity, and stress: 

* 2 * 3* 2 * /m \ 

p = ep , u a = e u a , v — e v , r = er. (21) 

Now we write the equations and boundary condition in the scaled vari- 
ables. In what follows, we drop the "*" superscript on all re-scaled variables. 
First, the continuity equation, V • u = 0, 

d d d 

-(h 2Ul ) + — (h lU2 ) + — (hhiv) = 0, (22) 



dx± dx2 dy 
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then, the Stokes momentum equation, Vp = V 2 u: 



ei dp e 2 dp 1 dp 1 
+ — h-e 3 - 



d (h, d 



h\dx\ h 2 dx 2 e dy h\h 2 \ dx\ \h\dx\ 
2 d (hi d 



+ 



^ d_ ( d_ 

dx 2 \ h 2 dx 2 ) dy \ 1 2 dy 



(tiie-i + u 2 e 2 + eve 3 ) . (23) 



The boundary conditions now take the following forms. 

• First, the no-slip boundary condition ([J) is 

Ui = , on y = . (24) 

• Now consider the zero-shear-stress boundary condition at y — rj. First, 
from (0), express the components of the normal unit vector n in terms 
of the scaled variables: 

cn a = -e(l - er)k a >)— , cn 3 = (1 - er)k x ) (1 - er)k 2 ) , (25) 

where c = \tx X t 2 \ is the constant of normalisation. Furthermore, 
from H p599], the components of the non-dimensional (symmetric) 
deviatoric stress tensor, r = (Vu + Vu T ), become, upon scaling: 



n 2 

r 33 



2e 



1 du c 
m n dx 



+ 



-u n > 



m\m 2 



+ O (e 2 ) 



1 du\ 1 du 2 mi 2 m 2 1 

H —Ui —u 2 



m 2 dx 2 mi dxi mim 2 
du a 
dy 

dy 



mim 2 



+ O (e 2 ) 



+ ek a u a + O (e 2 ) 



(26) 



Thus to order e, the only contributing terms in the boundary condition 
t • t a = is r a 3n 3 , namely 



dy 



+ ek a u a + O (e 2 ) = 



(27) 



In order to write down the normal stress boundary condition on the 
free surface, equation @, we first need to express the free-surface mean 
curvature R to order e: from ( [T7D the scaled free-surface mean curvature 
is 

k = K + 6 K 2 r] + eV 2 r? + O (e 2 ) , 
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where k = k± + k 2 and k 2 = kf + k\. Next, the normal surface stress 
component to order e 2 in scaled form is (using summation) 

-p + eniHjTij = -p + 2e — + O e . 

\ay mi ay m 2 ay ) v ' 

Hence we now write (||) as 

p = —k — eK 2 r] — eV 2 ?7 + O (e 2 ) , on y — rj ; (28) 

viscous stresses have no influence on the normal stress to this order in 
this scaling. 



3.2 Perturbation solution 

We now find a solution of these equations by assuming a perturbation ex- 
pansion of the unknown fields in terms of the small parameter e. We write 
each component in the following expansion: 



V = V 

(0) + ev (D + e 2 v (2) + . . . ^ (2Q) 

p = + ep {1) + e 2 p (2) H . 



Then, at the leading order, we find the following equations governing u_ 
v (o) an dp(°): 



(o) 

a i 



op (0) 

dy 

d 2 u^ 1 dp<® 
dy 2 m a dx a ' 



(30) 



-(7712^1 ) + q ( m l u 2 ) + "^1^2" 



with the boundary conditions 

u (o) = v (o) _ q at y = o , 



at 27 = r/, (31) 



= — k oX y = rj . 

The solution of these equations is readily found to be the expected locally 
parabolic flow driven by pressure gradients induced by curvature variations: 

P<°> = -K, 
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u^^-—K a (\y 2 -r ] y) , (32) 
(°) = V 2 k (|y 3 - \ m 2 ) - |Vk • V V y 2 . 



v 

At the next order of perturbation: 

dy 



0. 



cty a?/ m Q y ox a ox a ) 

d , (i). <9 , (i), dv {1) 



(m2u\ ) + — — {m\u 2 ) + m\m2- 



dx\ dx2 dy 

with the boundary conditions 

= = o a t y = , 

M 1 ) 



-^«i 0) at 2/ = ?;, (34) 



= — (« 2 f7 + V 2 7/) at y = rj 



We find 



— {k 2 T] + V 2 i]) ia (\y 2 -rjy)- -^-n ja y 3 (35) 



m a ' ^ ' 6m a 

lyS _ l_ vy 2 + I, 



/1Q 1 o 12\ 

«,a(g2/ ^ 2^ + 2^ 



3.3 Conservation of mass 

The expressions in the previous subsection show how the velocity and pres- 
sure fields respond to free-surface and substrate curvature. Such flow will 
thin the film in some places and thicken it in others. Conservation of fluid 
then leads us to an expression for the evolution of the film's thickness rj by the 
driven flow. Initially we work with non-dimensional but unsealed variables 
and coordinates before returning to scaled quantities. 

Consider a small volume above a patch of the substrate extending across 
the fluid layer from y = to y = r), see Figure ||. It is bounded by the 
substrate, the instantaneous free surface and the coordinate surfaces x%, x± + 
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dxii x 2 , x 2 + dx 2 . The rate at which fluid leaves this volume is expressed to 
first order of infinitesimal quantities dx\ and dx 2 as the sum of the terms 

J o [uih 2 dx 2 ] Xl+dxi dy - [uih 2 dx 2 } xi dy , (36) 

for the surfaces x\ and X\ + dx%, corresponding terms for the surfaces x 2 and 
x 2 + dx 2 , and the term 

/ u ■ ndS = \h\h 2 v — r/ \h 2 Ui — r/ 2 h\u 2 \ dx\dx 2 = hih 2 -^-dxidx 2 , (37) 
Jt l > j fit 

for the bounding free surface where we have used the kinematic boundary 
condition (jl5|) . Upon dividing by dx\dx 2 and in the limit of dx\dx 2 — > 0, we 
deduce 

dt] d d 

where the Q^s are the components of the total flux of fluid over a position 
on the substrate, defined as 

rv 

Q(t,xi,x 2 ) = Qiei + Q 2 e 2 = / ((1 - k 2 y)uie-i + (1 - hy)u 2 e 2 )dy . (39) 

Jo 

After division of (^) by mim 2 , recognize on the right-hand-side the surface 
divergence of the flux vector Q expressed in the orthogonal curvilinear coordi- 
nate system (xi, x 2 ) of S. This yields the divergence form of the conservation 
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of mass equation: 

(1 _ ^0(1 - Vh)^; = — f^-KQO + ^-( mi Q 2 )) = - V • Q. 

at m\m 2 \ax\ ox 2 J 

(40) 

We find an approximation of the flux vector Q by taking into account the 
thin layer characteristics of the flow over S as determined earlier. 

Returning to the rescaled quantities introduced in §§|3.1|, we now deter- 
mine the perturbation coefficients of the flux vector Q = e 2 Q(°) + e 3 Q^ 1 ' ) + • • • 
to be 

Q<°> = f rfV + ufe 2 )dy = ^ + 2 ) = \rf Vk , (41) 

jo \nii m 2 / 

Q (1) = jT ((u? - hyuf^e, + (4 1} - ^ 0) )e 2 ) dy 

= \rfV{K 2 r] + V 2 7]) + lr] 4 K- Vk-^kVk, (42) 

where k = ki + k 2 , k 2 = k 2 + k 2 , and K = k a 5 a p(e a : e^) is the curvature ten- 
sor (which would not be diagonal in a general orthogonal coordinate system 
of the substrate S). 

Finally, the corresponding evolution equation for 77, the so-called lubrica- 
tion approximation, is given by 



1 - er}ki)(l - er)k 2 )^- = -|e 3 V ■ h 3 V« - er/^nVn + |er/ 4 K ■ Vk 

+0 (e 5 ) , (43) 

where k = k + en 2 r] + eV 2 r/ + O (e 2 ). 

We express this equation in the more convenient form 



-±e 3 V- 

dt 3 



r/ 2 CV« - |er/ 4 (/d — K) • V«] + O (e 5 ) , (44) 



where ( = r/ — ^enrf + ^e 2 kik 2 r] 3 is proportional to the amount of fluid in the 
film lying "above" a patch of the substrate. Rewriting this in terms of un- 
sealed, non-dimensional quantities gives (0) when the free-surface curvature 
is approximated by an expression such as ([1] 



4 Centre manifold analysis: Gravity and in- 
ertia 

The previous analysis is adequate to give a basic model for the dynamics of 
the film thickness. However, it is not clear how to reasonably modify such an 
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analysis for generally curved substrates in the presence of more complicating 
physical effects such as gravity. In this section we turn to the more powerful 
centre manifold theory to extend the previous model in order to include a 
gravitational body force and to determine the leading effect of inertia. The 
systematic nature of the centre manifold approach lends itself well towards 
computer algebra and hence its results provide a valuable check on the hu- 
man algebra reported in the previous section. In future work, this approach 
will lead to the correct modelling of initial conditions and forcing, as has 



been done for dispersion in pipes and channels |L5|, e.g.], and perhaps small 



substrate roughness — see the review in for an introduction. This section 
also serves as an example of the analysis that may be undertaken if other 
physical effects are needed in the lubrication model of the dynamics. 

4.1 Basis of the centre manifold 

Centre manifold theory H assures us that the long-term dynamics near a 
fixed point of a dynamical system may be accurately modelled by a "low- 
dimensional" system which is based upon the linear picture of the dynamics, 
see the review |[26|| . Here the lubrication model ([I]) is of low-dimension in 
comparison to that of the incompressible Navier-Stokes equations (0-0). The 
first task is to establish the linear basis of the centre manifold. 



Analogous to the approach developed by Roberts in |23|, |24], we retain 
the small parameter e introduced in § |3.1| and scale the curvatures and met- 
ric coefficients according to (|18HT9|). By considering e negligible we, in the 



preliminary linear analysis, restrict attention to large-scale flows on a flat 
substrate. However, we do not explicitly scale the unknown fields, p and u 
as in (|2lD , because in this approach their scaling naturally arises from the 
dynamical equations during the course of the analysis rather than being im- 
posed at the outset — another virtue of this approach. Further, we treat the 



gravitational forcing as a small effect. By substituting, as in [23], the Bond 
number b = f3 2 , then adjoining the trivial equations 

de 83 

s = 0, and ^ = 0, (45) 



to the scaled versions of the fluid equations (|^-§D and their boundary con- 
ditions, we treat all terms involving products of u, e and 3 as perturbing 
nonlinear terms (as in the centre manifold analysis of the unfolding of bifur- 
cations |H). 

Such nonlinear terms are discarded in the preliminary linear analysis to 
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leave linear equations 

d\i 

V ■ u = , and 11— + Vp - V 2 u = , (46) 

at 

where, since e = 0, in effect these differential operators are the Cartesian 
operators appropriate to a flat substrate. These equations are to be solved 
with boundary conditions of u = on S (y = 0) and 

du dv 

~dy =P ~ 2 dy = ° ony = r l' ( 47 ) 

and the linearized kinematic condition 

drj 

— v — , on y = n . 

dt y ' 

All solutions of these linear equations are composed of the decaying lateral 
shear modes, v = p = 0, u a = c a sm(£7ry / (2r))) exp(At) for odd integers £, 
together with r\ = constant, u a = v = p = 0. In these modes the decay-rate 
in time is A = — £ 2 tt 2 /(4r] 2 7l), except for the last mentioned mode which, 
as a consequence of fluid conservation, has a decay-rate A = 0. So linearly, 
and in the absence of any lateral variations on a flat substrate, all the lateral 
shear modes decay exponentially quickly, on a time-scale of 7Zi] 2 , just leaving 
a film of constant thickness as the permanent mode. This spectrum, of 
all eigenvalues being strictly negative except for a few that are zero, is the 
classic spectrum for the application of centre manifold theory; the Existence 
Theorem 1 in || assures us that the nonlinear effects just perturb the linear 
picture of the dynamics so that in the long-term all solutions of the full 
nonlinear system are dominated by the slow dynamics induced by gravity 
and large-scale lateral variations in the film thickness and curvature. The 
Relevance Theorem 2 in assures us that these dynamics are exponentially 
attractive and so form a generic model of the long-term fluid dynamics of the 
film. With the caveat that strict theory has not yet been extended sufficiently 
to cover this particular application, the closest being that of Gallay || and 
Haragus [TT| (but also see jpTjl), we apply the centre manifold concepts and 
techniques to systematically develop a low-dimensional lubrication model of 
the dynamics of the film. 

Having identified the critical mode associated with the zero decay-rate, 
the subsequent analysis is straightforward. The usual approach is to write the 
fluid fields v(t) = (u\, u 2 , v, p), as a function of the critical mode rj (effectively 
equivalent to the "slaving" principle of synergetics Jl0|)- Instead of seeking 



asymptotic expansions in the "amplitude" of the critical mode [20, 15, 23 
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e.g.], we apply an algorithm to find the centre manifold and the evolution 
thereon which is based directly upon the Approximation Theorem 3 in || ^] 



and its variants, as explained in detail by Roberts in P5f . An outline of the 
procedure follows. 

We seek solutions for the fluid fields as 

v(t) = V(t?) , such that ^ = G(r], V(rj)) , (48) 

where dependence upon the constant parameters (e, (3) is implicit, and where 
G is the right-hand side of the re-scaled version of the kinematic condi- 
tion (|T5|). The aim is to find the functional V such that v(t) as described 
by ([48]) forms actual solutions of the scaled Navier-Stokes equations; this 
ensures fidelity between the model and the fluid dynamics. Suppose that 
at some stage in an iteration scheme we have an approximation, V and 
G = G(r],V), and then seek a correction, V, so that we obtain a more 
accurate solution to the governing equations. Substituting 

v(t) = V + V , such that ^ = G 

into the scaled Navier-Stokes equations then rearranging, dropping products 
of corrections, and using a zeroth-order approximation wherever factors mul- 
tiply corrections, we obtain a system of equations for the corrections which 
is of the form 

BV' = R, (49) 

where B is the linear operator seen on the left-hand side of (f4l)|-44~7|) , and, 
most importantly, R is the residual of the scaled Navier-Stokes equations 
using the current approximation. After solving this equation to find the 
correction V, the current approximation V and G is updated. The iteration 
is repeated until the residual of the governing equations, R, becomes zero to 
some order of error, whence the centre manifold model will be accurate to 
the same order of error (by the Approximation Theorem 3 in ||). 

A computer algebra program^ was written to perform all the necessary 
detailed algebra for this physical problem. A very important feature of this 
iteration scheme is that it is performed until the residuals of the actual gov- 
erning equations are zero, to some order of error. Thus the correctness of the 

2 The computer algebra package REDUCE was used because of its flexible "operator" 
facility. The source code is publicly available via the URL http : //www . sci . usq . edu . au/- 
"robertsa or by contacting the second listed author. At the time of writing, information 
about reduce was available from Anthony C. Hearn, RAND, Santa Monica, CA 90407- 
2138, USA. E-mail: reduce@rand.org 
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results that we present here is based only upon the correct evaluation of the 
residuals and upon sufficient iterations to drive these to zero. Thus the key 
to the correctness of the results produced by the computer program is the 
proper coding of the fluid dynamical equations in the chosen curvilinear coor- 
dinate system. Upon obtaining the code, these can be seen in the computed 
residuals within the iterative loop. Also note that because the thickness of 
the film is continuously varying in space and time, it is convenient to work 
with equations in terms of a scaled vertical coordinate Y = y/ry(x, t) so that 
the free surface of the film is always Y = 1. We also let X a = x a , T = t com- 
plete the new coordinate system. Then space and time derivatives transform 
according to 

d Id d d Y di] d d d Y Or/ d 
dy = r]dY' &~ ~ dX^ ~ ~q^X^8Y 1 di ~ df ~ ^dfdY ' 

However, the fluid equations are not rewritten in this new coordinate system 
because the computer handles all the necessary details of the transformation. 



4.2 The general lubrication model 



Based upon scalings introduced in §3.1, we run the reduce computer al- 
gebra program that computes velocity and pressure fields for the flow and 
the evolution equation for the film's thickness. We find that the evolution 
equation may be written in the coordinate-free form| 



dt 



rfCVk - \eq A (Kl - K) • Vk 



€7] A (kI + §K) ■ g s + eg y rfVrj 



0(e 5 + b< 



(50) 



where g s and g y are respectively the components of the gravitational unit vec- 
tor tangent and normal to the substrate. Note that in the special orthogonal 
curvilinear system used in the derivation of the model, 



kI-K 



h o 
ki 



and kI + |K 







2/ 

fei + %k 2 



lh + k 2 



One example of the application of the model ( j50|) is to the dynamics of 
thin films on vertical cylinders or fibres as explored by Kalliadasis & Chang 



3 (e p + b q ) is used to denote terms, z, for which z/(e p + b q ) is bounded as (e, b) — » 0. 
The upshot is that z — e m b n is O (e p + b q ) iff m/p + n/q> 1. For example, an expression 
accurate to O (e 5 + b 2 ) retains all terms of the form e m b n for 2m + 5n < 10. 
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12 



They find that a thin film may form into drops on the fibre, or may 
saturate into solitary waves depending upon the film thickness. The model 
they use, derived by Frenkel |7] , is approximately the version of ( pUD specific 
to axisymmetric flow around a cylinder of radius a, namely 



dt 



ld_ 

3dz 



9V 



+ v 2 C 



l 



drj d 3 r] s 



(a + 77) 2 dz dz 3 







(51) 



in our nondimensionalisation, where z measures axial distance and here C = 
n+i] 2 /2a. The differences between our systematic approach and theirs is that 
their model, (2) in 



12 



does not conserve fluid, while our model (pip does 
and also has some higher order corrections to account better for the curvature 
of the substrate. Further, (|50|) specialised to a cylindrical substrate will also 
descibe any non-axisymmetric dynamics of interest. 

Returning to the general model fl5T|), it is derived solely under the as- 
sumptions that curvature and free-surface slopes, as measured by e, and the 
gravitational forcing, measured by the Bond number b, are perturbing influ- 
ences. The application of centre manifold theory places no restrictions upon 
their relative magnitudes — we do not have to insist on b ~ e or any other 
such relation between these two independent parameters. Provided there are 
no "run away" instabilities, the model is valid over any scaling regime where 
both parameters are small. In particular, the model is valid as the tangential 
and normal gravitational forcing vary widely around a curving substrate. 

Also note that the model (|50|) was derived without placing any overt 
restriction upon the Reynolds number, it was treated as an O (1) constant. 
That the model actually turns out to have no Reynolds number dependence 
just confirms the Stokes flow nature of these lubrication dynamics. Higher- 
order analysis, as also seen in |23|, §4], shows that inertia first appears at 
O (e 6 + b 3 ) for fluid films. Thus inertia is formally negligible. However, if the 
Reynolds number is large enough to be significant, then we have determined 
that the following O (Tie 6 + Tib 3 ) terms should be included in the evolution 
equation (|50|) : 



dt 



K 

T 

2Kb 



15 



V- — (V// • Vk) Vk - ?rV (Vk • Vk) + — (V' 2 k)Vk 

3 7 63 

V • 7] 6 (g fl • Vr] Vk + g s Vr] ■ Vk) - -L- V (g s • Vk) 



+ -^g s V 2 K + -^gyKVK 



21 
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2Kb 2 



15 



■V ■ 



g s -V77) + ^(13Kl-9K).g s ^ 



(52) 



Note that for a flat substrate {k\ — &2 = 0), these Reynolds number correc- 
tion terms reduce to those derived by Benney and Atherton & Homsy . 
We suggest that these terms are more likely to be used to estimate the error 
in the lubrication model ([50]) when applied to moderate Reynolds number 
flows. Thus these expressions may be used to indicate when a more sophisti- 
cated dynamical model, such as a two-mode model || |23| , is required in order 
to resolve the inertial instabilities of fluid films at higher Reynolds numbers. 
Also recall from the previous subsection that the time-scale of decay of the 
lateral shear modes is lZrj 2 and so the other criterion for validity of the model 
is that the time-scales exhibited in a simulation are significantly longer than 
this. 

Of course, an order of magnitude argument will give a global estimate, 
over space and time, of the influence of inertia. However, one may make 
a perfectly satisfactory a priori order of magnitude assessment, but if in a 
simulation an instability grows, then (^) will detect if it grows too large. 
Thus (|2|) supplements order of magnitude estimates by giving focussed local 
estimates of the errors in an actual simulation. 



5 Example film flows 

In this section we use the preceding lubrication models to demonstrate some 
of the flow effects caused by substrate curvature. 

5.1 Qualitative effect of substrate curvature 

Equation (|44]) shows that, to leading order of perturbation, the flow is 
driven by substrate curvature gradients, curvature caused by film thickness 
variations is generally smaller unless the substrate is comparatively gently 
curved. To understand the effect of curvature qualitatively, consider a two- 
dimensional flow on a one-dimensional substrate with given curvature k. De- 
note x = Xi as substrate arclength so that the scale factor m x = 1 and there 
are no variations in X2- Then to leading order and in the absence of body 
forces, the flow is governed by the first-order partial differential equation for 

7]: 

(l-^ = -||-(^). (53) 
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This has a characteristic solution 

i = j ?7 = ~^ K xx~. • (54) 

1 — K7] a 1 — KT) 

• Wherever the substrate curvature k has a local minimum, say at x = xq, 
then 

k x < (x < x ) , n x > (x > x ) , and k xx > . 

Thus according to the characteristic equations the flow is thinned in 
the neighborhood of x = Xq. Indeed, since rj oc —n xx r] 3 , the film thick- 
ness r] typically decreases as t^ 1 ^ 2 at a point of minimum curvature 
(k xx > 0). As shown schematically in Figure^, this thinning applies in 
the neighbourhood of minimum absolute curvature for interior coating 
flows, and of maximum absolute curvature for exterior coating flows. 

• Conversely, the film thickens in the neighbourhood of a point of max- 
imum curvature as also shown in Figure |^. The characteristic solu- 
tion ( (54| ) predicts a film thickness that blows up in finite time. The 
role of higher-order terms in the model equations is to smooth out such 
singularity. 
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5.2 Corner flow 

Consider here the two-dimensional fluid dynamics examined by Schwartz & 



Weidner [g9| consisting of the flow which thins a film around an outside 
corner as shown in Figure |j. The fluid is taken to have surface tension 
a = 30dynes/cm, viscosity \i = lpoise, and density p = lgm/cm 3 . Initially 
the film is a uniform 0.01cm thick, corresponding to a Reynolds number of 
1Z = 0.003, around the corner which has a radius of 0.1cm. As in the previous 
subsection, the most convenient way to parameterise the one-dimensional 
substrate is in terms of the arclength x. In this case our dynamical model ([!]) 
reduces to 



d(_ f dr] l d 

dt~ [l KV) dt~ *dx 



2 dk 
dx 



O (e 5 ) , (55) 



where k = k + k 2 t] + and k = ki(x) is the curvature of the substrate. In 



dx 

contrast, the model of Schwartz & Weidner (the SaW model) is 



dr] y d 
dt 3 dx 



d ( d 2 r) s 



V 3 ir I k 



dx \ dx 



2 



(56) 



The differences between the models are that ours includes more terms in the 
curvature. In particular, ours conserves fluid whereas the SaW model does 
not. Indeed, in the course of the numerical simulations of the film flow shown 
in Figure |], the SaW model lost about 2% of the fluid whereas ours lost none 
to computational error. For thicker films the difference is more marked. 

A numerical scheme to simulate a fluid film via these equations is straight- 
forward. We approximated both (|55"D and ( |56l) by finite differences on a spa- 
tial grid with N = 97 points and, because the dynamics are stiff, we employed 
an implicit integration scheme with time step At = 0.00115s. The numer- 
ical scheme uses second-order centred differences in space and time but, in 
the interests of speed, the nonlinear coefficients are only computed from the 
earlier time. By varying the size of the space-time grid we determined that 
these parameters give a numerically accurate simulation. 

Shown in Figures |5| and § are comparisons between the predictions of 
the SaW model and ours during the simulation of the thinning of the film 
around the corner. Observe that the SaW model and ours are quantitatively 
different: the SaW model predicts a more rapid thinning of the film around 
the corner, and a slower thickening of the film away from the corner. For 
example, the thickness at x = and t = 3s for our model is only reached by 
the SaW model at time t ~ 300s. For quantitative accuracy the thin film 
flows need the extra curvature terms employed in our model 
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Figure 4: The initial film of fluid (dot-dashed) of constant thickness around 
a corner (solid) evolves over a long time to thin around the corner (dashed). 
(Unlike Schwartz & Weidner |29|, Fig. 7], the thickness of the film has not 
been exaggerated.) 
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Figure 5: film thickness for the flow around the corner shown in Figure |4] at 
time t = 3s: dashed is our model (|55|) ; dot-dashed is the SaW model (|56|) . 
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Figure 6: film thickness for the flow around the corner shown in Figure |] at 
time t = 300s: dashed is our model (|55|) ; dot-dashed is the SaW model (|56D. 
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Strictly, this corner flow is unsuitable for the lubrication approximation 
because the sharp change in curvature from the flat wall to the round corner 
generates arbitrarily large gradients of the curvature k. These gradients are 
compounded at higher order. For example, a computation of the Reynolds 
number terms in the first line of (|52]) for this problem reveals a relative 
magnitude, when compared to drj/dt, of typically 10007?. ~ 3. Thus such 
high-order effects may be significant near the sharp change in curvature. In 
practise, perhaps we should view the sharp change as an element of roughness 
in an otherwise smooth transition between the wall and the round corner. 
Alternatively, we could view it as a discontinuity between two smooth (con- 
stant) curvature domains and thereby seek matching conditions via similar 
arguments to those employed in ET[ for boundary conditions. 



5.3 Flow around a torus 

Here we discuss the evolution of the flow over the surface of a torus with 
tube of radius R 2 and centreline of radius R\. We use the following parame- 
terisation: 

X\ = (Ri + R2 cos 9) cos 4> , 

X 2 = (R t + i? 2 cost?) sin 0, (57) 
X 3 = R 2 sin9. 

Take R 2 < R\ to avoid self-intersection, and denote the coordinates x\ = (p 
and x 2 = 9 with < <p < 2tt and < 9 < 2tt. Then this parameterisation 
generates an orthogonal curvilinear coordinate system with corresponding 
orthonormal basis vectors: 

e i = e 4> = — sin 0i + cos 0j , 

e 2 = e e = — sin ^(cos 0i + sin 0j) + cos 9k , (58) 
e 3 = ei x e 2 = cos #(cos cf>i + sin 0j) + sin 9k , 

in terms of the unit vectors of the conventional Cartesian coordinate system, 
and with the corresponding surface metrics 

mi = m<p = Ri + R 2 cos 9 , m 2 = m e = R 2 . (59) 

The chosen coordinate system also generates lines of principal curvature as 
the parameter curves. We then obtain the following principal curvatures and 
mean curvature: 

cos# 1 1 R 1 + 2R 2 cos9 

K\ = , Ko = , K= . oO 

R 1 + R 2 cos9 , R 2 ' R 2 R 1 + R 2 cos9 K ' 

The mean curvature of the substrate is maximum at the inner rim of the 
torus (9 = n), and minimum at the outer rim (9 = 0). Hence we expect the 
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fluid layer to thicken around the inner rim, and to thin around the outer, 
solely due to surface tension effects. 

We numerically solve the model (fH|). For simplicity, we seek axisymmet- 
ric solutions, that is solutions independent of the angle <fi around the rim of 
the torus, and so the film thickness depends only upon 9, the angle around 
the tube, and t. We try the following form for r) 

N-l 

V (t,6)= £ a n {t) cos(n0), (61) 

n=0 

which guarantees the periodicity of the solution and imposes symmetry across 
the plane z — 0. The ordinary differential equations for the coefficients a n 
are found by a Galerkin method where they are determined by making the 
corresponding residual error orthogonal to the N basis functions, cos(n#), in 
the usual L 2 norm. As a check, we confirm that the total volume enclosed be- 
tween the free surface of the fluid and the toroidal substrate remains constant 
in time. Two numerical simulations were done. 

• First, we start with an initially uniform layer of thickness r] = 0.1 on 
a torus with Ri = 2 and R% — 1- The corresponding flow towards the 
inner rim of the torus is shown by the evolution of the film thickness 
in Figure [7J 

• Second, we simulated the flow evolving from a strip of fluid placed 
around the outer rim of the torus. Again, as shown in Figure ||, the 
fluid flows around to the inner rim. 

On the torus, the importance of inertia effects on this lubrication flow are 
estimated by the ratio between typical values of the right-hand side of 
and that of fl52|): dominantly 

■ , n ■ ■ Kt] 4 
inertia : lubrication m 2 . 

As may be expected, a thicker film or a more sharply curved torus are more 
likely to be affected by such higher-order influences. Note that for flows 
inside a toroidal tube, the thinning of the liquid layer occurs around the 
inner rim (6 = tt). Three-dimensional simulations of coating flows exhibiting 
transversal instabilities (along the ^-direction) similar to those found in [39] 
will be published elsewhere. 

Acknowledgements: RVR would like to acknowledge the many helpful 
discussions with LW. Schwartz. His work was completed while spending a 
sabbatical leave at USQ whose warm hospitality is also gratefully acknowl- 
edged. 
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Figure 7: evolution of flow on the surface of a torus with N = 15 terms in 
Galerkin approximation. Ri = 2, R 2 = 1, i]q = 0.1: (a) rj(t,9) shown at 
t = 0, 100, 300, 1000, 3000, 10000; (b) initial, intermediate (t = 1000) and 
late (t = 10000) stages of the film on a cross-section of the torus. 



Roy, Roberts & Simpson, February 9, 2008 



5: Example film flows 



30 




Roy, Roberts & Simpson, February 9, 2008 



§5: Example film flows 



31 



References 

[1] R.W. Atherton and G.M. Homsy. On the derivation of evolution equa- 
tions for interfacial waves. Chem. Eng. Comm., 2:57-77, 1976. 

[2] G.K. Batchelor. An Introduction to Fluid Dynamics. CUP, 1979. 

[3] D.J. Benney, Long waves on liquid films. J. Math & Physics, 45:150-155, 
1966. 

[4] R.A. Cairncross, L.F. Francis and L.E. Scriven. Predicting drying in 
coatings that react and gel — Drying regime maps. AICHE Journal, 
42:55-67, 1996. 

[5] J. Carr. Applications Of Centre Manifold Theory, volume 35 of Applied 
Math Sci. Springer- Verlag, 1981. 

[6] H.C. Chang. Wave evolution on a falling film. Annu. Rev. Fluid Mech., 
26:103-136, 1994. 

[7] A.L. Frenkel. Nonlinear theory of strongly undulating thin films flowing 
down vertical cylinders. Europhys Lett, 18:583-588, 1992. 

[8] Th. Gallay. A center-stable manifold theorem for differential equations 
in Banach spaces. Commun. Math. Phys, 152:249-268, 1993. 

[9] J.B. Grotberg. Pulmonary flow and transport phenomena. Annu. Rev. 
Fluid Mech., 26:529-571, 1994. 

[10] H. Haken. Synergetics, An introduction. Springer, Berlin, 1983. 

[11] M. Haragus. Model equations for water waves in the presence of surface 
tension, preprint, 1995. 

[12] S. Kalliadasis and H.-C. Chang. Drop formation during coating of ver- 
tical fibres. J Fluid Mech, 261:135-168, 1994. 

[13] V.G. Levich. Physicochemical Hydrodynamics. Prentice-Hall, 1962. 

[14] C.C. Mei. The Applied Dynamics of Ocean Surface Waves, World Sci- 
entific, 1989. 

[15] G.N. Mercer and A.J. Roberts. A centre manifold description of con- 
taminant dispersion in channels with varying flow properties. SIAM J. 
Appl. Math., 50:1547-1565, 1990. 



Roy, Roberts & Simpson, February 9, 2008 



§5: Example film flows 



32 



[16] J. A. Moriarty, L.W. Schwartz, and E.O. Tuck. Unsteady spreading of 
thin liquid films with small surface tension. Phys. Fluids A, 3:733-742, 
1991. 

[17] J. A. Moriarty and L.W. Schwartz. Effective slip in numerical calcula- 
tions of moving-contact line problems. J. Eng. Math., 26:81-86, 1992. 

[18] J. A. Moriarty and L.W. Schwartz. Dynamic considerations in the clos- 
ing and opening of holes in thin liquid films. J. Colloid Interface Sci., 
161:335-342, 1993. 

[19] P. Morse and H. Feshbach. Methods of Mathematical Physics. McGraw- 
Hill, New York, 1953. 

[20] A.J. Roberts. The application of centre manifold theory to the evolution 
of systems which vary slowly in space. J. Austral. Math. Soc. B, 29:480- 
500, 1988. 

[21] A.J. Roberts. Boundary conditions for approximate differential equa- 
tions. J. Austral. Math. Soc. B, 34:54-80, 1992. 

[22] A.J. Roberts. The invariant manifold of beam deformations. Part 1: the 
simple circular rod. J. Elasticity, 30:1-54, 1993. 

[23] A.J. Roberts. Low-dimensional models of thin film fluid dynamics. Phys. 
Letts. A, 212:63-72, 1996. 

[24] A.J. Roberts. An accurate model of thin 2D fluid flows with inertia on 
curved surfaces. In P. A. Tyvand, editor, Free- surface flows with viscosity, 
Advances in Fluid Mechanics Series. Comput Mech Pub, 1996. 

[25] A.J. Roberts. Low-dimensional modelling of dynamics via computer 
algebra. Comput Phys Comm, 100:215-230, 1997. 

[26] A.J. Roberts. Low-dimensional modelling of dynamical systems, 
preprint, USQ, Feb 1997. 

[27] G.J. Roskes. Three-dimensional long waves on a liquid film. Phys. 
Fluids, 13:1440-1445, 1969. 

[28] K.J. Ruschak. Coating flows. Annu. Rev. Fluid Mech., 17:65-89, 1985. 

[29] L.W. Schwartz and D.E. Weidner. Modeling of coating flows on curved 
surfaces. J. Engrg. Math., 29:91-103, 1995. 



Roy, Roberts & Simpson, February 9, 2008 



§5: Example film flows 33 

[30] L.W. Schwartz, D.E. Weidner, and R.R. Eley. An analysis of the effect 
of surfactant on the leveling behavior of a thin liquid coating layer. 
Langmuir, 11:3690-3693, 1995. 

[31] L.W. Schwartz, R.A. Cairncross and D.E. Weidner. Anomalous behavior 
during leveling of thin coating layers with surfactants. Phys. Fluids A, 
8:1693-1695, 1996. 

[32] J.B. Sweeney, T. Davis, L.E. Scriven, and J. A. Zasadzinski. Equilib- 
rium thin films on rough surfaces. 1. Capillary and disjoining effects. 
Langmuir 9:1551-1555, 1993. 

[33] J.J. Stoker. Differential geometry. Wiley-Interscience, New York, 1969. 

[34] S.P. Timoshenko and S. Woinowsky-Krieger. Theory of Shells and 
Plates. McGraw Hill, 1959. 

[35] E.O. Tuck and L.W. Schwartz. A numerical and asymptotic study of 
some third-order ordinary differential equations relevant to draining and 
coating flows. SI AM Review, 32:453-469, 1990. 

[36] B.W. van de Fliert, P.D. Howell, and J.R. Ockendon. Pressure-driven 
flow of a thin viscous sheet. J. Fluid Mech., 292:359-376, 1995. 

[37] D.E. Weidner and L.W. Schwartz. Contact-line motion of shear-thinning 
liquids. Phys. Fluids, 6:3535-3538, 1994. 

[38] D.E. Weidner, L.W. Schwartz and R.R. Eley. Role of surface tension 
gradients in correcting coating defects in corners. J. Colloid Interface 
Sci., 179:66-75, 1996. 

[39] D.E. Weidner, L.W. Schwartz and M.H. Eres. Simulation of coating 
layer evolution and drop formation on horizontal cylinders. J. Colloid 
Interface Sci., 187:243-258, 1997. 



Roy, Roberts & Simpson, February 9, 2008 



